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The controllability of a network is a theoretical problem of relevance in a variety of 
contexts ranging from financial markets to the brain. Until now, network controllability 
has been characterized only on isolated networks, while the vast majority of complex 
systems are formed by multilayer networks. Here we build a theoretical framework 
for the linear controllability of multilayer networks by mapping the problem into a 
combinatorial matching problem. We found that correlating the external signals in 
the different layers can significantly reduce the multiplex network robustness to node 
removal, as it can be seen in conjunction with a hybrid phase transition occurring 
in interacting Poisson networks. Moreover we observe that multilayer networks can 
stabilize the fully controllable multiplex network configuration that can be stable also 
when the full controllability of the single network is not stable. 


Most of the real networks are not isolated but interact with each other forming multilayer structures [U |2] . For 
example, banks are linked to each other by different types of contracts and relationships, gene regulation in the cell 
is mediated by the different types of interactions between different kinds of molecules, brain data are described by 
multilayer brain networks. Studying the controllability properties of these networks is important for assessing the 
risk of a financial crash [3111], for drug discovery [5| and for characterizing brain dynamics [BHTQ]. Therefore the 
controllability of multilayer networks is a problem of fundamental importance for a large variety of applications. 

Recently, linear dlHISI and non-linear [20H30j approaches are providing new scenarios for the characterization of 
the controllability of single complex networks. In particular, in the seminal paper by Liu et al. [12] the structural 
controllability of complex networks has been addressed by mapping this problem into a Maximum Matching Problem 
that can be studied using statistical mechanics techniques [5T1136) . Other works approach the related problem of 
network observability or target control (SSj which focuses on controlling just a subset of the nodes. Despite the 
significant interest in network controllability, all linear and non-linear approaches for the controllability of networks 
are still restricted to single networks while it has been recently found that the multiplexity of networks can have 
profound effects on the dynamical processes taking place on them [551144j . For example, percolation processes that 
usually present continuous phase transitions on single networks can become discontinuous on such structures 
and are characterized by large avalanches of disruption events. 

Here, we consider the elegant framework of structural controllability m and investigate how the multilayer structure 
of networks can affect their controllability. We focus on multiplex networks, which are multilayer networks in which 
the same set of nodes are connected by different types of interactions. Multiplex network controllability is studied 
under the assumption that input nodes are the same in all network layers, thus mimicking the situation in which 
input nodes can send different signals in the different layers of the multiplex but the position of the external signals 
in the layers is correlated. 

We show that controlling the dynamics of multiplex networks is more costly than controlling single layers taken 
in isolation. Moreover, the controllability of multiplex networks displays unexpected new phenomena. In fact these 
networks can become extremely sensible to damage in conjunction with a discontinuous phase transition characterized 
by a jump in the number of input points (driver nodes). A careful investigation of this phase transition reveals that 
this is a hybrid phase transition with a square root singularity, therefore in the same universality class of the emergence 
of the mutually connected component in multiplex networks [DIMIIIII- The number of driver nodes in the multiplex 
network is in general higher than the number of driver nodes in the single layers taken in isolation. Nevertheless 
the degree correlations between low-degree nodes in the different layers can affect the controllability of the multiplex 
network and modulate the number of its driver nodes. Moreover, a fully controllable configuration can be stable in a 
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The replica nodes are either all driver 
nodes or none is 


FIG. 1: Control of a multiplex network. The controllability of a duplex network (multiplex with M — 2 layers) can be 
mapped to a Maximum Matching Problem in which the unmatched nodes (indicated with a white circle) are the driver nodes 
of the duplex network. Here we have indicated with red thick links the matched links and by black thin links the unmatched 
links. 


multilayer network even if it is not stable in the isolated networks that form the multilayer structure. 


RESULTS 

We consider multiplex networks [T] in which every node i = 1,2,..., N has a replica node (i, a) in each layer a and 
every layer is formed by a directed network between the corresponding replica nodes. We assume that each replica 
node (i, a) is characterized by a different dynamical variable xf S R and that each layer is characterized by a possibly 
different dynamical process. We consider for simplicity a duplex, i.e a multiplex formed by two layers {A, B} where 
each layer a € {A, B} is a directed network. The state of the network at time t is governed by a linear dynamical 
system 


dX{t) 

dt 


gx + icu. 


( 1 ) 


in which the 2Wdimensional vector X(t) describes the dynamical state of each replica node, i.e. Xi = xf and 
Xjv+i = xf ior i = 1,2,..., N. The matrix is a 2N x 2N (asymmetric) matrix and Ail is a 2N x P matrix. They 
have the following block structure 


G 


9 ^ 0 
0 


1C = 


Q \ 

0 ^ 


( 2 ) 


in which are the N x N matrices describing the directed weighted interactions within the layers and K°‘ are the 
N X P“ matrices describing the coupling between the nodes of each layer a and < N external signals. The 
latter are represented by a vector u(t) of elements u~^ and 7 = 1, 2... P = P^ + P^. Here we consider the concept 
of structural controllability mnn] that guarantees the controllability of a networks for any choice of the non-zeros 
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entries of Q and /C, except for a variety of zero Lebesgue measure in the parameter space. Therefore each layer of the 
duplex networks can be structurally controlled by identifying a minimum number of driver nodes, that are controlled 
nodes which do not share input vertices. If different replicas of the same node can be independently controlled, then 
the controllability properties of the multiplex network factorize and each layer can be studied as if was taken in 
isolation [iiiiiiiTniiin]. Liu et al. [12] showed that in a single layer the minimum set of driver nodes can be found 
by mapping the problem into a matching problem. In real multiplex networks however nodes are usually univocally 
defined and share common properties across different layers, therefore we make the assumption that each node of the 
duplex network is either a driver node in each layer or it is not a driver node in any layer. The problem of finding 
the driver nodes of the duplex network can be thus mapped into a maximum matching problem in which every node 
has at most one matched incoming link and at most one matched outgoing link, with the constraint that two replica 
nodes either have no matched incoming links on each layer or have one matched incoming link in each layer (see 
Figure]^. This problem can be studied, using statistical mechanics techniques, such as the cavity method and the 
Belief Propagation (BP) algorithm. Following [T2[T2], we consider the variables s“- = 1,0 indicating respectively if 
the directed link from node (i, a) to node (j, a) in layer a = A, B \s matched or not. In order to have a matching in 
each layer of the duplex the following constraints have always to be satisfied 

E <1. E < 1- (3) 

iedZj 

where d^i is the set of replica nodes (j, a) in layer a that are reached by directed links from (i, a) and dfj is the set 
of replica nodes {i,a) in layer a that point to (j, a). In addition, we impose that the driver nodes in the two layers 
(the unmatched nodes) are replica nodes, i.e. 


E =5 = E 4 - (■>) 

iEd^j iGdfj 

In this formalism, computing the maximum matching corresponds to minimize an energy function B = No = Nno 
where No is the number of unmatched replica nodes associated to each matching. The energy E for a given matching, 
can be expressed in terms of the variables Sij as 


E = E •s')- w 

a j \ i^dZj j 

In order to study this novel statistical mechanics problem, we derived the BP equations [S31 ESI ES] (see Methods and 
Supplementary Material) valid in the locally tree-like approximation, as described for the case of a single network 

problem in [II m EH El EH- 


DISCUSSION 

The controllability of multiplex networks displays a rich phenomenology, coming from the interplay between the 
dynamical and the structural properties of multiplex networks. Here we characterize the controllability of multiplex 
networks with different degree distribution and with tunable level of structural correlations. 

Phase transition in Poisson duplex networks- We consider duplex networks in which the two layers are realizations 
of uncorrelated directed random graphs characterized by Poisson distributions for in-degree and out-degree with same 
average value c, i.e. (fc^) = {k^^f) = (fc^) = {kgut) = c- In Figure]^ we report the average rescaled number of driver 
nodes n_D as function of the average degree c computed from the solutions of Eqs.Q on single instances and from the 
graph ensemble analysis. A comparison with two independent layers with the same topological properties shows that 
the controllability of a duplex network is in general more demanding in terms of number of driver nodes than the 
controllability of independent single layers, in particular for low average degrees. In addition, a discontinuity in the 
number of driver nodes at c = c* = 3.2223 ... marks a change in the controllability properties of duplex networks that 
is not observed in uncoupled networks. This is due to a structural change in the solution of the matching problem, 
in which a finite density of zero valued cavity fields emerges. A careful investigation (see Supplementary Material ) 
reveals that this is a hybrid phase transition with a square root singularity, therefore in the same universality class of 
the emergence of the mutually connected component in multiplex networks [HEIIIIT]. 

In correspondence to this phase transition the network responds non trivially to perturbations. This is observed by 
performing a numerical calculation of the robustness of the networks. Following HU we classify the nodes into three 


4 



O' ^ ^ ' . 

0 1 2 3 4 5 6 



c 


c 


FIG. 2: Controllability of Poisson duplex networks with average degrees = 

c. In panel A the fraction no of driver nodes in a Poisson duplex network with = c, 

plotted as a function of the average degree c. The points indicate the average BP results obtained over 5 single realizations of 
the Poisson duplex networks with average degree c and N = 10^, the solid line is the theoretical expectation (the error bar, 
indicating the interval of one standard deviation from the mean, is always smaller or comparable to marker size). The dashed 
line represents twice the density of driver nodes for a single Poisson network with the same average degree. 

In panel B the densities ric, rir and no respectively of critical redundant and ordinary nodes are shown as functions of c for the 
same type of duplex networks with N = 10^, where each point is the average over 100 different instances. 

In both panels the dot-dashed vertical line indicates the phase transition average degree c* = 3.222326106.... 


categories: critical nodes, redundant nodes and ordinary nodes. When a critical node is removed from the (multiplex) 
network, controllability is sustained at the cost of increasing the number of driver nodes. If the number of driver nodes 
decreases or is unchanged, the removed nodes are classified as redundant and ordinary respectively. Figure shows 
that the fraction of critical nodes reaches a maximum at the transition, revealing an increased fragility of the duplex 
network to random damage with respect to single layers. While an abrupt change in the number of driver nodes can 
result from a small change in the network topology, it is important to stress that the non-monotonic behavior of these 
quantities around the critical average degree value could be interpreted as a precursor of the discontinuity. 

In a duplex network formed by directed Poisson random graphs with different average degree in the two layers (i.e. 
(kZi) = {Kut) = Ca) a similar discontinuous phase transition is observed (see Figure]^. Nevertheless we checked 
that this discontinuous phase transition is not occurring for every multiplex network structure (see Supplementary 
Material). 

Effect of degree correlations on the controllability of duplex networks - We consider a model of duplex network in 
which the replica nodes of the directed random graphs in the two layers have correlated degrees. In particular, we 
consider a case in which only the low in-degree nodes (nodes with in-degree equal to 0,1, 2) are correlated (replica 
nodes in different layers have same degree with probability p) and a case in which the in-degrees of the replica nodes 
are equal with probability p independently of their value (see Supplementary Material for details). The controllability 
of the network is affected by these correlations as shown in Figure]^ In fact, the number of driver nodes no decreases 
as the level of correlation increases. In duplex networks with Poisson degree distribution, low-degree correlations 
modify both the position of the hybrid transition and the size of the discontinuity. Once the replica nodes with low 
in-degree are correlated, a further correlation of the remaining replica nodes does not substantially change the number 
of driver nodes. This result confirms that structural controllability is essentially determined by the control of low 
degree nodes [l6] . 

Stability of the fully controllable solution - A fully controllable solution, in which a single driver node is necessary to 
control the whole duplex network, exists if the minimum in-degree and the minimum out-degree are both greater than 
I in both layers. This solution of the cavity equations gives the correct solution to the maximum matching problem 
describing the controllability of multilayer networks only if no instabilities take place. The stability conditions are then 
found by imposing that the Jacobian of the systems of equations derived by the cavity method has all its eigenvalues 
Ai of modulus less than one, i.e. |Ai| < I. In random duplex networks with the same degree distribution in the two 
layers, the fully controllable solution is stable (see Supplementary Material for the details of the derivation) if and 
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FIG. 3: Phase diagram of the controllability for a Poisson duplex networks with average degrees = 

^kA,out\^ _ and = cb- The color code indicates the density of driver nodes no = E/N in the multiplex 

network. 




FIG. 4: The effect of the degree correlation between replica nodes in different layers on the controllability 
of multiplex networks. Correlations between the low in-degrees (Id) and correlations between any in-degree node (td), 
parametrized by p, affect the fraction of driver nodes in the network no, both in the case of Poisson networks with the same in 
and out average degree c across the two layers (Panel A) and in the case of scale-free networks with the same in and out degree 
distribution across the layers, given by P{k) oc k~'^ and minimum in/out degree 1 (Panel B). When p = Q there is no degree 
correlation between replica nodes in different layers. The BP data are shown for networks with N = 10"*, and are averaged 5 
times for panel A and 20 times for panel B. 
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FIG. 5: The fully controllable solution can be stable for the multiplex network also if it is not stable for the 
single layers taken in isolation. Histogram of the number of networks that out of 100 realizations have No driver nodes. 
The results obtained for the control of a duplex networks and its two layers are compared. The duplex networks are formed by 
two scale-free networks with N = 10^ and = pB,in ^ ^ 2 and oc k ' for k > 2, with 

7 = 2.3, the networks have minimum in-degree equal to 2 and minimum out-degree equal to 3 and P^’®"(2) = P^’®"(2) = 0.3. 


only if 


pout ^ 2 ) < 


2 -!))„’ 


( 6 ) 


for a = A, B. On single networks it was instead recently found m that the fully controllable configuration is only 
stable for 


P“(2) < 


(k) 


2 (fc(fc-l)), 


pout (^ 2 ) < 


(k) 


2{kik-l))^ 


(7) 


This implies that for multiplex networks with asymmetric in-degree and out-degree distributions it might occur that 
the fully controllable solution is stable in the multiplex network but unstable in the single networks taken in isolation 
(see Figure for the characterization of the controllability of a similar type of multiplex networks). Therefore a 
multiplex structure can help to stabilize the fully controllable solution. 

In conclusion, within the framework of structural controllability, we have considered the controllability properties 
of multiplex networks in which the nodes are either driver nodes in all the layers or they are not driver nodes in any 
layer. Our results show that controlling multiplex networks is more demanding, in terms of number of driver nodes, 
than controlling networks composed of a single layer. In random duplex networks with Poisson degree distribution, it 
is possible to observe a hybrid phase transition with a discontinuity in the number of driver nodes as a function of the 
average degree, that is phenomenologically similar to the emergence of mutually connected components. Close to this 
phase transition the duplex network exhibits an increased fragility to random damage. The existence of correlations 
between the degrees of replica nodes in different layers, in particular between low-degree nodes, has the effect of 
reducing the number of driver nodes necessary to control duplex networks. Finally, multiplex structure of networks 
can stabilize the fully controllable solution also if this solution is not stable in the single layers that form the multiplex 
network. 


METHODS 

The BP equations - The BP equations of this problem are derived using the cavity method [3311351136j as described 
for the case of a single network problem in [niiaisiiisiisi. The same approximation methods can be applied here. 
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as long as the structure of the interconnected layers is locally tree-like both within the layers and across them. Under 
the decorrelation (replica-symmetric) assumption, the cavity fields (or messages) and defined on the 

directed links between neighboring nodes (i,a) and {j,a) in the same layer a = A, B satisfy the zero-temperature 
limit of the BP equations, also known as Max-Sum equations. 




— max 


— 1, max 


h 


cx. 

k^i 




— max 


max hf.i , — max hu.j 
kGd^i\j kaa^j 




— max 


max h ?.,, — max ht.i 
ked^i\j ^ kea^j 


( 8 ) 


in which the fields are defined to take values in the discrete set {—1, 0,1} and here and in the following we use the 
convention that the maximum over a null set is equal to —1 (see Supplementary Material for details). In terms of 
these fields, the energy E in Eq. becomes 


N 


i V 

E = — max 

— 1, max 

/cG i 


Q J a 1 r a 

O' i—1 


“ (i,i) 



i=l,N 


max 


0 , max ht 
kea^i 


max h 
kGd^i 


B 

'k—>-i 


(9) 


BP equations over ensemble of networks- Let us consider the case of uncorrelated duplex networks in which the 
degree of the same node in different layers are uncorrelated and there is no overlap of the links. In each layer a = A, B 
we consider a maximally random network with in-degree distribution P“’*"(A:) and out-degree distribution P“’°“*(/c). 
At the ensemble level, each link of (the infinitely large) random network forming layer a has the same statistical 
properties, that we describe through distributions Va{h°‘) and Va{h°‘) of cavity fields that are defined on the support 
of Eqs. 1^ i.e. 


Vo.(h°‘) = wf5{h°^-l)+wf5{h°^+ l)Pwf5{h°‘), 

V^{h°‘) = w‘f5{h°‘-l)+w^5{h°‘+ l)+wf5{h°‘), (10) 


where a = A, B and where the probabilities wf,W 2 , are normalized wf wtf = I as well as the probabilities 
wf ,11)2,11)^ that satisfy the equation wf wlf = 1. The cavity method at the network ensemble level is also 

known as density evolution method [36j . 


It is useful to introduce the generating functions Gq’™^°“*(z), and G' 


a,inf out 


G. 


oc^infout 
0 


(z) of the multiplex network as 


W = 


/out^k)z’^, g'“>™/°“*( 2 ) = with a = A, B. In this way, we can 


derive recursive equations for the probabilities {ref }i^i^ 2,3 and {ref } 2 ^i, 2 , 3 , that are the analogous of the BP equa¬ 
tions for an ensemble of uncorrelated duplex networks 


Wn = 


Wo = 


wf 


Gf’“(zc^) [l-G^’™(l-u;f) 

1 - Gu (I - ) -h Gu (1 - )Go ’ (^2 ) 

G B.in/ B\ 1 

1 ) [1 - Go (1 - ) 

I ^B,in 


-I B\ I y^rsan/-, B\r^ 

1 - ’ (1 - rci ) -h ’ [1- Wi )G, 


A,in 
0 


W) 


( 11 ) 


with Wg = 1 — wf — W 2 , and = I — wf — wtf. The energy E and the entropy density s of the matching problem 
can be also expressed in terms of the {w“}i=i_ 2,3 and {ref }i=i, 2.3 (see Supplementary Material for details). 

Hybrid transition for Poisson duplex network- Here we consider the case of two Poisson networks with the same 
in/out average degree. In other words, we consider the situation in which = c. 






















We notice that the BP equations can be rewritten to form a closed subsystem of equations for wi and W 2 (see 
Supplementary Material for details), 


xci = hi{wi,W 2 ) = e 


— r.p. 


1 — e 




W 2 = h 2 {wi,W 2 ) = 1-e +e e 


(12a) 

(12b) 


from the solution of which the remaining quantities can be determined. 

The value c* of the average degree c at which the discontinuity in the number of driver nodes nu observed is found 
by imposing that Eqs. (12) are satisfied together with the condition 


|J|=0, 


(13) 


with J indicating the Jacobian of the system of equations (12|. Imposing that Eqs. (121 and condition (13) are 


simultaneously satisfied, the solution c* = 3.222326106... is found. Eor c < c* we observe that W 3 = ws = 0. At c* 
we observe a discon tinuit y i n bo th W 3 and W 3 , but for c > c* the functions hi{wi,W 2 ) and h 2 {wi,W 2 ) are analytic, 
and analyzing Eqs. (12a I — (126) we obtain the behaviour of the order parameters W 3 and W 3 for c> c* 


ics — u >3 oc (c— 
W 3 — w^ oc (c— 


(14) 


showing that the transition is hybrid. 
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Supplemental Material ” Control of Multilayer Networks” 

A. Introduction 

This Supplemental Material is structured as follows. 

In Sec. II we define the problem of structural controllability of multiplex networks, focusing on the case of a duplex 
network. Moreover we define the driver nodes, as the set of nodes that, if stimulated by an external signal, can drive 
the dynamical state of the network to any desired configuration. 

In Sec. Ill we map the problem of structural controllability of a duplex network to a Maximum Matching Problem, 
and we derive the Belief Propagation (BP) equations determining the driver nodes, and their zero-temperature limit 
known as Max-Sum equations. 

In Sec IV we consider the controllability of uncorrelated duplex networks, characterizing the BP equations valid for this 
problem, the stability conditions for the solutions of the BP equations, and the entropy of the solutions. Moreover, we 
consider duplex networks formed by two Poisson layers and we characterize their hybrid phase transition. Finally we 
consider the controllability of duplex networks formed by layers with power law in-degree and out-degree distributions. 
In Sec. V we consider ensembles of duplex networks in which the in-degrees of replica nodes are correlated and we 
derive the BP equations assuming either that only the low in-degrees of replica nodes are correlated or that all the 
in-degrees of replica nodes are correlated. 


B. The structural controllability of a multiplex network 


We consider a multiplex network in which every node i = 1, 2,..., has a replica node in each layer and every layer 
is formed by a directed networks between the corresponding replica nodes [T] . We assume that each replica node can 
have a different dynamical state and can send different signals in the different networks (each layers is characterized 
by a different dynamical process). In this case the controllability of the multiplex network can be treated by control 
theory methods used for the single layers taken in isolation [ni m [13 HO]. Nevertheless here we will consider an 
additional constraint on the number of driver nodes. In fact we impose that corresponding replica nodes are either 
driver nodes in all layers or they are not driver nodes in any layer. 

We consider for simplicity a duplex, i.e a multiplex formed by two layers where each layer is formed by a directed 
network. We call the two layers layer a = A, B. We consider a linear dynamical system determining the network 
dynamics 

^=0X(t)+/Cu(t), (15) 

in which the vector X(t) describes the dynamical state of each replica node in the duplex, and has 2N elements. The 
first set of N elements represents the dynamical state xf of node i in layer A (i.e. Xi = xf for i = l,2...,iV), 
while the elements Ajv+i represent the dynamical state of the node i in layer B, and are given by X^+i = xf for 
j = 1, 2,..., iV. The matrix ^ is a 2N x 2N (asymmetric) matrix and the matrix /C is a 2N x M matrix. The matrices 
Q and JC have the following block structure 


g = 




x = 


Q \ 

0 )' 


( 16 ) 


where g°^ with a = A, B are the N x N matrices describing the directed weighted interactions within each of the 
networks in the two layers and Ar“ are the N x matrices describing the interaction between the nodes of the 
network a and the M“ < N external signals for layer a. The external signals are indicated by the vector u{t) of 
elements and 7 = 1,2 ... M = . 


Given block structure of both matrix Q and X described in Eq. (16 1 , the problem of duplex network controllability 
defined by Eq. (15), can be exactly recast into the problem of control ability of the single layers that form the duplex 
network. 

Here we adopt the framework of structural controllability [n] aimed at characterizing if a given duplex network is 
controllable when the non-zero matrix elements of g and X given by Eq. (16 1 are free parameters. A duplex networks 
in which the linear dynamics described by the Eqs. (151 and (16) take place, is structurally controllable if both layers 
a = A,B are structurally controllable. 

Each layer a is structurally controllable if for any choice of the free parameters in g°‘ and AT", except for a variety of 
zero Lebesgue measure in the parameter space, the Kalman’s condition is fulfilled m- Since structural controllability 
only distinguishes between zero and non-zero entries of the matrices 5 “ and K°‘ , a given directed network in layer a is 
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structurally controllable if it is possible to determine the input nodes (i.e. the position of the non-zero entries of the 
matrix K°‘) in a way to control the dynamics described by any realization of the matrix g°‘ with the same non-zero 
elements, except for atypical realizations of zero measure. In practice, a single network can be structurally controlled 
by identifying a minimum number of driver nodes, that are controlled nodes which do not share input vertices in both 
layers. In their seminal paper |12j . Liu and coworkers showed that on single networks this control theoretic problem 
can be reduced to a well-known optimization problem: their Minimum Input Theorem states that the minimum set 
of driver nodes that guarantees the full structural controllability of a network is the set of unmatched nodes in a 
maximum matching of the same directed network. Their result for a single network remains valid for the duplex 
network described by Eqs. (151 — (16). Therefore the structural controllability of duplex networks, in the absence 
of further constraints can be mapped to a Maximum Matching problem defined on the single layers of the duplex 
networks. Here nevertheless, we consider a further constraint to be imposed on the driver nodes, which enforce a 
new type of dependence between the layers of the duplex. In particular we impose that the replica nodes (i, a) with 
a = A, B and a given index i, are either both driver nodes or neither is a driver node. This implies that these two 
replica nodes are either both linked to independent and external signals or none of them is connected to external 
signals. 


C. The Maximum Matching Problem for the Controllability of Duplex Networks 

1. Mapping duplex controllability into a constrained Maximum Matching Problem 


In order to build an algorithm able to find the driver nodes of a duplex network we consider the variables s“- = 1, 0 
indicating respectively if the directed link from node (z, a) to node (j, a) in layer a = A, B is matched or not. In the 
two layers of the duplex network we want to have a matching, i.e. the following constraints must always be satisfied 
for a = A,B, 


E »3 

j&dp 

< 1, 

(17a) 

H *3 

jGd°i 

< 1, 

(17b) 


where here and in the following we indicate with 9“ the set of nodes j that are pointed by node i in layer a and with 
the set of nodes j pointing to node i in layer a. In addition we impose that the driver nodes in the two networks are 
replica nodes, i.e. in the matching problem either two replica nodes are both matched or both unmatched. Therefore 
the variable s“- satisfy the following additional constraints 


E 4 = E 


(18) 


iGd^j 




Finally we need to minimize the number of driver nodes in the multiplex network. Therefore we minimize the energy 
E of the problem given by 





with 


^r = i- E 




jGdp 


The energy E is given by the number Nd of driver replica nodes in the duplex network by 

E = = Nti d . 


(19) 


( 20 ) 


( 21 ) 


2. Derivation of the BP equations at finite inverse temperature fi 

We consider here the Maximum Matching Problem defined in Sec | C 1| The goal is to find the configuration of the 
variables } associated to every directed edge z —>■ j in layer a, such that the energy E given by the number of 
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driver replica nodes in the duplex network is minimized provided that the conditions given by Eqs. (Hii, @ are 
satisfied. Introducing as an auxiliary variable the “inverse temperature” /3 we cast this problem into a statistical 
mechanics problem where our first aim is finding the distribution P({sij}), parametrized by the inverse temperature 
(3, and given by 




a-PE 


N 


n n 


E Ml- E 




j^d+i 


jed-i 


E E 

iGd^j i^d^j 


4 


( 22 ) 


where 9{x) = 1 for a: > 0 and 9{x) = 0 for a; < 0 , S{x) is the Kronecker delta, and where Z is the normalization 
constant, that corresponds to the partition function of the statistical mechanics problem. Subsequently, we plan to 
perform the limit /? —>■ oo in order to characterize the optimal (i.e. the maximum-sized) matching in the network 
satisfying Eqs. ( [l8| . The free-energy of the problem F(/3) is defined as 

/3F(/?) = -lnZ, (23) 


and the energy E is therefore given by 


E = 


d[PFm 
dp ■ 


(24) 


The distribution T’({sij}) on a locally tree-like network can be (approximately) estimated by the cavity method in 
the replica symmetric assumption (i.e. by deriving Belief Propagation equations) [T^ [T51 [5TH5B] . In this respect, in 
each layer a of the duplex network, we define two probability marginals on each directed link, one going in the same 
direction of the link P°^j{sij) and one in the opposite direction P°^j{sji). The BP equations for these quantities are 






1 


V: 


E 


1 


1- E 


s^k exp 


-Hi- E 


^ik 




n pA^Atk) 

ked^i\j 


(25a) 


E 


T)^ . 


1- E 


^ki 


E 

{s^^}\k^d^i 


E 


A 


k£d^i 


E E ' 

iGd^j i^d^j 


X n n 

k^d-^i\j k^d^i 


(25b) 


P^ = 

\^JlJ 


1 


- 




Ski 


E 




{s^J|/cG9d* 


1- E 


^ki 


kGd^i 


E E ^ 


X Y[ pAAsZ) n PAiAt) 

k£d^i\j k^d^i 


(25c) 


where and are normalization constants. The probability marginals {P^j{sfj), P^j{s'j^)} can be 

parametrized by the cavity fields hf_^, and defined by 




i—^j\ ij / l+exp[^/3/i; 




(26) 


In terms of the cavity fields (or messages), Eqs. (25) reduce to the following set of finite temperature BP equations, 

I 


= -^log|e- 


E - 

kGd^i\j 


PK 


■ 


= -jlog 


1 


= -^log 


1 


E 

k^d^i\j 

E ' 

kedt^i\j 


C.PK 


oPht 


(27a) 

(27b) 

(27c) 
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The free energy F and the energy E = of the model are given respectively by 


N 


-/^^ = EE 


ol i—1 


In e-« + 




and by 


N 


E 


= EE 


a. i—1 


0-0 _ 




O 0 K 




P.^^k- 


E 1 M 1 + E 


o 0 hi 


i=l,N 


E E 111 (1 + , 






kGd^i k'Gd^i 






-E 


E. 


kGd^i 


o 0 ht^ 




1 + J2ked^i Efe' 


GOfi 


,- 6 ' 


0 hZ_ 


E E 

Ct: <i,j>a 






■i+h° 


1 + e 


0 ihz 


. +/l9 


(28) 


(29) 


3. BP Equations for /3 —^ 00 


The BP equations in the limit ^ 00 are derived from the Eqs. (27). In the limit /? —)■ 00 the solution is expressed 


in terms of the helds hf_o^j or hf_^j sent from a node {i,a) to the linked node {j,a) in layer a = A, B. The cavity 
fields have a simple interpretation as messages between neighboring replica nodes m ■ = hf_^j = 1 means 

“match me”, hf_^j = hf_^j = — 1 means “do not match me”, and = 0 means “do what you want”. The 

zero-temperature BP (or Max-Sum) equations determining the values of these fields in the limit /3 —>■ 00 are are given 
by 


hf_^j = — max 




max hu.i, — max h?., 
kGd^i\j kaa^j 


max h?.,, — max hf., 
kea^i\j ked^j 


in which the fields are defined to take values in the discrete set {1,0, —1} and we defined the maximum over a 
set equal to —1. It follows that for fcf= 0 we have = —1 and for = 0 we have hf_^j = —I. 

The energy E can also be expressed in terms of these fields and is given by 


(30c) 
null 


N 


iV 

^ = -EE 

a i—1 

+EE 

a <i,j> 


max 


max 


— 1 , max /lEi 

fcG9“i 

0 . A-h°‘ 


— max 
i=l,N 


0, max /lEj -I- max h?., 
k^d^i k^d^i 


(31) 


where < i,j > indicates pair of nodes that are nearest neighbors in the network and where we take the maximum 
over a null set equal to -1. 
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D. Controllability of uncorrelated multiplex networks with given in-degree and out-degree distribution 

1. Cavity equations for an uncorrelated multiplex network ensemble 

Let us consider the case of uncorrelated duplex networks in which the degree of the same node in different layers 
are uncorrelated and there is no overlap of the links. In each layer a = A, B we consider a maximally random 
network with in-degree distribution and out-degree distribution At the ensemble level, each link 

of (the infinitely large) random network forming layer a has the same statistical propert ies, that we describe through 
distributions Va{h‘^) and of cavity fields that are defined on the support of Eqs 


30 


i.e. 


= w‘fSih°‘-l)+w^6ih°‘+ l)+w^d{h°‘), 


(32) 


where a = A^B and where the probabilities w“, , iCg are normalized wf = 1 as well as the probabilities 

Wi,W 2 tW^ that satisfy the equation -|- -I- wf = 1. The cavity method at the network ensemble level is also 

known as density evolution method [36] . 

It is useful to introduce the generating functions Gq’*"''^°“*(z), and Qf multiplex network as 

G“’“(z) = ^P“’“(fc)z^ 


1 V ; 


{k)z‘ 

h 

ik)z' 


k-l 


G“’°“‘( 0 ) = yjT^p' 


(33) 


with a = A,B. In this way, we can derive recursive equations for the probabilities {w“}i=i, 2,3 and {thf }i=i, 2 , 3 , that 
are the analogous of Eqs. 30 for an ensemble of uncorrelated duplex networks 

wf = G“’°“‘(u;2 “), 

= [l-Gr“‘(l-^ir)], 

W 3 = I-W 3 - W 2 , 

W 3 = I-W 3 - W 2 , 

= Gf’“(u;^) [l-G^’“(l-u;f) 


1 - Gi ’ (1 - ’ (1 - )Go ’ (^^2 ) 


A,in / 


''tS B\ 1 /^/i,272/1 A\ 

Wl = Gl ) 1 - Go ’ (1 - ) 


u>f = 


1 ^BAn/^ B\ I /^B.in/^ B\ryA,in ( A\ 

1-Gg (1-u;i) + Gg (1-u;i)Go’ (u>2 ) 


(34) 


The energy E of the matching problem can be also expressed in terms of the {ref }i=i, 2,3 and (ref }i=i, 2,3 giving 


E = Y {Go“’°“‘ «) - [1 - Go“’°“*(l - ref)] }-{[!- G^’“(l - ref)][l - Go^’“(re 2 ^)] 

OL 

+ [1 - G^’™(1 - ref )][1 - Gf’“(re^)]| + ^ (fc“)i„ [ref(1 - ref) -b ref (1 - ref)]. 


(35) 


.S. Stability condition 


The Eqsj^ might have multiple solutions. In order to evaluate the stability of these solutions, using a method 
already used in the context of single networks [Hiin] here we compute the Jacobian of the system of Eqs. p4|) and 
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impose that all its eigenvalues have modulus less than one. We avoid to consider and because they influence 
only the number of null eigenvalues (4 eigenvalues upon 12). The 12 x 12 Jacobian matrix becomes 8 x 8 and it can 
be decomposed in four 4x4 blocks 


with 


iJii = 


J = 


Hll i/l .2 
H 2.1 H 22 


VG^’“( 1 -O(l-Go^>“(u;f)) 


0 

0 G^’°“'(l - ref) 

Gf’“(ref)(l-Go^’“"(l-ref)) 0 

0 0 


0 Gf’°“*(ref) \ 

A,out/-I 


^A,out 

0 
0 
0 


, (36) 


H 2 A = 


0 0 0 0 \ 

0 0 00 ' 

Gf’“(zef)(fc)^,„Gf’“(l-0 0 0 0 

0 G^”(l-0(A:)^,„Gf’“K) 0 0/ 


(37) 


and 


H, 2 = 


'1^2= Gf’“K)(fc)3,„Gf’“(l-0 


0 0\ 
0 0 
0 0 


Gf’“(l-0(A:)^,,„Gf’“(ref) 0 Oj 


B,in/ B\ 


H2,2 = 


0 0 G 2 ’ (w^) \ 

0 Gf’°“*(l-wf) 


\G 2 "’“(l- 0 (l-Gf’“K)) 


Gr’“(ref)(l-Gf'“(l- 0 ) 

0 


,out \ 

0 
0 
0 


/ 


(38) 


(39) 


Here the generating functions ( 7 “’®"/°“* ^re given by Eqs. (33) and the generating functions 

and G 2 ’°“*(a;) are defined as 


Gr”(^) = 

^ /tn 

Gr^\z) = 

, )out 


(40) 


Of particular interest is the characterization of the stability of the solution wf = wf = wf = ^2 = 0 and = 1, 

corresponding to the full controllability of the network, a conhguration with E = Njj = 0. This solution emerges for 
-^^"(l) = -^^*(1) = 0 for a = A,B. Therefore if the minimum in-degree and the minimum out-degree are both greater 
than one, the analysis at the ensemble level is consistent with the full controllability of the network. Nevertheless this 
solution might be not stable. By analyzing the Jacobian J for uif = = tef = '*^2 = 0 and = 1, we can 

determine under which condition the full controllability solution is stable. The Jacobian matrix, in this case simplify 
significantly and is given by 


( 0 


0 

0 

W)E~ 

0 


\ 


0 

0 

0 


0 0 

0 0 

0 0 

0 0 

0 0 

0 0 

0 0 


2P°“‘(2) 

0 

0 

0 

0 

0 

0 

0 


0 

0 

0 

0 

0 

0 

0 


0 0 

0 0 

0 0 

0 0 

0 0 

n 

0 0 

0 0 


0 ^ 

0 

0 

0 

2pr*(2) 

{k^)out 


0 

0 

0 


/ 


( 41 ) 
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Four eigenvalues of J are zero, the other four have degenerate modulus, therefore the stability conditions are 

{kA)out 

(fc®)m {k^)out 


(42) 


When P™{k) = P^*{k) = P^{k) = Pg^*{k) = P{k) we have just one stability criterion solution and it reads 


P(2)< 


{k? 


2{k{k-l)) 

We observe here that on the single layers a = A,B the full controllability solution, is instead only stable [16] for 


(43) 


^ (fc<-(fc°-l)),„Pr(2) 

(/c“)„ {k^)out 

^ (fc<-(fc°-!))„„, J^»(2) 

{k°^)out {k°^)in 


(44) 


It follows that for duplex networks in which both layers have the same in-degree and out-degree distributions, i.e. 
P“’*”(fc) = P“’°“‘(A:) the stability of the full controllability solution on single layers is the same as the stability on 
the duplex network. Nevertheless, for duplex networks formed by layers in which the in-degree distribution and the 
out-degree distribution are not the same there can be cases in which for the duplex network the fully controllable 
solution is stable while for the single layers it is not stable (See main text for discussion of this phenomenon and 
simulation results). 


3. Entropy 


In order to evaluate the number of maximum matchings, here we evaluate the entropy of the ground state solutions 
in the case of uncorrelated layers. The entropy density is given by Sg = Sq/N and can be computed by expanding the 
free energy at low temperatures /(/3 —)■ oo) = eg — sg//3 -I- 0(l//3^). This involves the study of the evanescent parts 
of the cavity field. Therefore we assume that the field can be written as 


7 _ 1 I In Uc. 

ILq, — i- ^ ^ 

h = -1 - 1 - Aie 

IL(^ X “T ^ 

h - 


for the peak around h = \ 
for the peak around h = —1 
for the peak around h = 0 


^0 = 1 +^ 
^ _1 ^ ^ 


i ^ InTo 

'^Oc ^ 


0 


for the peak around 
for the peak around 
for the peak around 


h = 1 
h = -l 
h = 0 


From the BP equations, and the equations for P{ha) and P{ha) we can obtain the relation between the probability 
distributions .A“(i^ct), (/x^), ^3 (yc), and the distributions ^“(i>q,), (Mo ):-^3 (7a); given by 
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A^y) = 


AU^i) = 


t + 1 ) / [ri«i 5 (f.f)l < 

fc=o J Li=i J V ^ + T,i=l^^f 

E ^o(™ + i)( 7 )«)''(i-»?r-‘ 

k=l ^2 ) V^/ 


^3(7) = 


X I n<ixMf(x.") 

E E i +1) (7) (%“)‘(<)“-‘ 


X / n‘<^.”xi;(7f) < 7-=^ 

x' U=i J \ 2^i=i 7i 




00 / J -| \ 00 00 ^ ^ \ 

= ^=0 kB = lmB=kB 7^ / 

f Y[dnfA^{fif)Y[diyj^Af{i^l^) xs(i)'^ - ^ ^ 


_ 1 _L ,,A 

' Z-^i — 1 r^i 


4^(A) = ^Gf’“(l-w(^)^i2®(A:^)(wf)''"’ / Jl (/xf) 


1 / I 1 \ / -Ax 

^2 kA = lmA=kA Vfc y 


X f UdixfAiiiyf) dU^- I \ 

'’ [i=i J V Hi^xK) 


^E E 

2 m^—k^ 


/2A(^x4 ^ ^ Pg{k^){wif 

7 4 JUB 


X / n ) n Tif (Aif) 5 (A^ — - —pi—7 

*=i i=i V + Xli=i 

Z^i = l H'i 
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-1 / B\ 

^^(7) = ^ ^ piik^)rA{wi)\wiy 

^3 fcB=lmB=fcB ^ ^ 


n dlFAi (7f) '^ ( 7^ - ^ 


OO OO / /4 I 1 \ / A\ 

[m +^) -cjA/^A , \/ A\k( A\m^-k^ 

+i)U.jK)K) 




OO OO / 

^ ^ [i^l 2-1 


X (5 7 — 


1 _ 1 I I 

V 7 ^ ^ ^ 2=1 /2 y 

-j OO OO / /4 I 1 \ / A\ 

;7(i_Gf"(I-<)) i: i: 77/^7".-+i)f“)(»7‘Kr'-‘' 


= l rn^ — k^ 


X / M7^- fc. ^ 

i=i V z;i=i7* 


The free energy density /(/3) = F{/3)/N = eo — ^ + (!l(l//3^) with 


where sqj are given by 


So — So,o,A + So,a,B + So,6 + So,c,^ + So,c,B 




+ J^^outikXw^fln 1 + J2 p 


OO OO / \ ^ 

+ EE (7)(^3)"(^2r-"x“.*Min ^7. 

k—lm—k \2=1 



18 


So.b — 


C 30 oo / 

^ ^ i=l 


oo oo / k^ 


fc^ = l m^—k^ 

oo oo 


oo oo / 

(i-G^’“(i-^«f)) E E (7A)KEKr"”'^"^n(m-")inE7.^ 

oo oo / Q\ k^ 

—Im^ —i —1 
oo oo / 

j:pS,('‘‘‘)(^>?)‘’'E E (’rJK)‘"(i-K)r”“‘>™(">'*) 




'In l + E^^'^E 


K 


i=l i=l 


_ oo oo / B\ 

- k^ k^ — lm^—k^ 


/ h.B uA \ 

'Infl + E-fE^f 


i=l 2=1 
A 


oo oo / 


oo oo / k^ 

E E (7B)K)'^"Kr"-'"^-(n^^)ln l + E^^^E^f 

fcB = lm®=fc^ ^ V i=l i=l 


(53) 


So,c.Q = —(fc“)m {w“(w“ + W3 )ln j>“ + t(;“(w“ + Wg )lnz/“ 

+ «)“w2 ln(l + ln(l + J^/t) 

+ Z()“w3 ln7 + ^“zDg ln7 + Wg lig ln(l + 77)|. ( 54 ) 


4- Phase transition in the controllability of Poisson duplex networks. 

Here we consider the case of two Poisson networks with the same in/out average degree. In other words, we consider 
the situation in which = (fc®’”®®*) = c. The fraction no of nodes that are driver nodes 

of this duplex, is always larger than the double of the fraction of driver nodes in each of the layers taken in isolation 
(see Fig. Moreover we observe that there is a phase transition in the controllability of these duplex networks, 
indicated by a discontinuity of no for c = c* = 3.22233... (see Fig. [^. In order to derive these results, we assumed 
wf = wf for i = 1, 2,3 and wf = wf for i = 1, 2, 3. Therefore the zero-temperature BP equations at the ensemble 
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FIG. 6: Density of driver nodes £ = no for a duplex network composed by two Poisson networks with = 

(fc®'*") = = c is indicated with a solid red line and it clearly shows a phase transition for c = 3.22233.... In the 

dashed blue line we display the double of the number of driver nodes e = 2n_D for a single Poisson network with the same 

average degree c, indicating the fraction of driver nodes necessary to control separately the two layers. 


level (34) read. 


Wi 
W 2 = 

W3 = 
W3 = 
Wi = 

W 2 = 


— g-c(l-'lll 2 ) 


1 — Wi — W2, 

1 — Wi — W2, 

g-c(l-u,2) _ g-cu,ij ^ 

1 _ 


The energy E is given in this case by 


E = 2 


o-c(i-'*2) _ 1 _(_ 


- 2[1 - 


^][1 


_ „-c(l-1«2 


+ 2 c[tai(l — W2) + wi(l — ^2)] ■ 


We notice that the equations for wi and W 2 can be rewritten to form a closed subsystem of equations, 

—c(l —'IW 2 ) 


— c.e. 


Wi = hi{wi,W2) = e 
W2 = h2{wi,W2) = 1-e 


1 — e 

-c(l-™ 2 ) 


p-c(l-™ 2 ) 


(55) 

(56) 

(57a) 

(57b) 
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FIG. 7: Plots of the functions wi = h\{w\,'W 2 ) and W 2 = h 2 {wi,W 2 ) given by Eqs. (57al — (576l. The solution of the system 
of these two equations, corresponds to a crossing of the two curves. We show the emergence of two new solutions of this system 
of equations for c > c* = 3.22233 .... The critical point c* characterize an hybrid phase transition in the controllability of the 
duplex network. 





FIG. 8: Values of the probabilities {wi}i=i, 2,3 and 'u)i=i, 2,3 plotted as a function of the average degree c, for a duplex network 
formed by two Poisson layers with = c. These probabilities are calculated directly 

from BP results obtained over 5 single realizations these multiplex networks with average degree c and N — 10"^. 
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from the solution of which the remaining quantities can be determined. 

The value c* of the average degree c at which the discontinuity in the number of driver nodes nu observed in Figj^ 
occurs can be found by imposing that the two curves rDi = hi(wi,W 2 ) and W 2 = h 2 (wi,W 2 ) of the plane Wi,W 2 for 
c = c* are tangent to each other at their interception. These functions are plotted in Figure]^ where it is possible to 
observe that for c > c* the curves cross in three points while for c < c* they cross in one point, and at c = c* they 
are tangent to each other. The critical point c* is found by imposing that the Eqs. ( [57| ) are satisfied together with 
the condition 


|J|=0, 

with J indicating the Jacobian of the system of equations wi = hi{wi,W 2 ) and W 2 = h 2 {wi,W 2 ) given by 


(58) 


J = 


1 - 


dhi{wi,W2) 

dw\ 

dh2{wi,W2) 

dwi 


1 - 


dhi{wi,W2) 

dw 2 

dh2{'Wi,W2) 

dw2 


Imposing that Eqs. (57) and condition (58) are simultaneously satisfied, the solution c* = 3.222326106... is found. 
For c < c* we observe that = W 3 = 0. At c* we observe a discontinuity in both W 3 and W 3 , but for c > c* the 
functions hi{wi,W 2 ) and h 2 {wi,W 2 ) are analytic, and analyzing Eqs. (57a) — (576) we obtain the behavior of the 
order parameters W 3 and 163 for c > c* 

W 3 — oc (c— 

W 3 -W 3 oc (c-c*)^/^, (59) 

showing that the transition is hybrid. 


We further characterize this phase transition evaluating the number of maximum matchings, i.e. the entropy value 
of the ground state solutions in the case of two poisson uncorrelated layers. The entropy density s follows from Eq. 
(51) and it is plotted as a function of the average degree c in Eig. The entropy density presents a small jump at 
c* = 3.22233.. marking a change in the properties of the solutions. 


Here we want to modify the degree distribution of the duplex network characterized in this section, by changing the 
probability of nodes of low degree (degree 0 , 1 , 2 ) that have been shown to be essential to determine the controllability 
of single layers [TB]. Therefore we consider a duplex networks with degree distributions P™(k) = P^*{k) = P'^{k) = 
P^^{k) = P{k) and with minimum degree is 2. In particular we consider P{k) given by 


P[k) = 


0 

P(2) 


for k <2 

for k = 2 

for k G [3, 00 ] 


(60) 


with K indicating a normalization constant. In Fig. 10 (on the left) we show the phase diagram of this duplex network 
described by the dependence of the fraction of driver nodes ud on c and P{2). The dark grey area defines the region 
where the zero-energy solution is stable, hence in which to control a duplex one needs only an infinitesimal fraction 
of driver nodes, i.e. no = 0. These results are compared with the situation in which the two layers are controlled 


separately shown in Fig. 10 (on the right). The fraction of driver replica nodes of the duplex network is always larger 


that the double of the fraction of driver nodes in any single layer taken in isolation. Moreover the region in which the 
fully controllable solution is stable is the same for the duplex network, and for the single networks in the layers of the 
duplex network taken in isolation. This result is consistent with the theoretical expectations obtained in Sec. | D 2| In 
fact the in- and out-degree distributions of the two layers are the same. 


5. Controllability of scale-free duplex networks 


Following Sec. | D 4| we consider now the case of two uncorrelated layers composed by two power-law networks with 
p(^A.zn) _ p(^A,outj _ P(/c^’*”) = P(fc®’°“‘) = P{k) oc k~'^ and minimal degree m = 1. Similarly to the poisson 
case, the fraction tid of driver nodes of this duplex, is always larger than the double of the fraction of driver nodes in 
each of the layers taken in isolation (see Fig. 11). 

Moreover, low-degree nodes significantly affect the controllability of duplex networks formed by scale-free networks. 
We consider a duplex network with P“(fc) = P^*{k) = Pf^{k) = P'^*{k) = P{k) and P{k) given by 


[0 if fc = 1 

P{k) = I P{2) if k = 2 

I Kfc-T' if fc e [3, M] 


( 61 ) 
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FIG. 9: Entropy density s for a duplex network composed by two Poisson networks with = 

= c. At c* = 3.22233.., average degree corresponding to the hybrid phase transition, the entropy density s displays a 
finite jump. 


with K indicating the normalization sum and 7 > 2. We consider uncorrelated networks, therefore the cutoff M on 
the degrees of the nodes will be given by 


M = min(yiV, {[1 - P(l) - P(2)]fV}^/^'^”^^). 


(62) 


In other words, the cutoff M is given by the minimum between the structural cutoff of the network and the natural 
cutoff of the degree distribution. In Fig. 12 (on the left) we present the phase diagram of a duplex network displaying 
the fraction of driver nodes nu as a function of the parameters 7 and P{2). The dark grey area is associated 
with the stable zero-energy solution while outside this region, the minimum fraction of driver nodes necessary for a 
full duplex control follows the colorcode. We compare these results with the situation in which the two layers are 
controlled separately (on the right). We observe that the number of driver replica nodes in the duplex is always 
greater than the total number of driver nodes of the single layer taken in isolation, provided that the duplex network 
is not fully controllable. We note that for the degree distribution considered in this case, consistently with the 
theoretical results obtained in Sec. D2 we observe that the region for the stability of the full controllability solution 


for the duplex network is the same of the region for the stability of the full controllability solution in the single 
layers. Finally , in Fig. we compare our theoretical results for the ensemble of d uplex networks with degree 
distributions P^{k) = P^^^k) = P^{k) = P^^{k) = P{k) and P{k) given by Eq. (61), with those obtained by the 
message-passing (BP) algorithm, finding a good agreement (Eq. [4^ returns a limit value for P(2) equal to 0.181947). 
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FIG. 10: On the left: the density of driver nodes n_D as a function of the para met ers c and P{2) is plotted for duplex networks 


with PTik) = Pa (k) = PM'ik) = Pb (k) = P{k) and P{k) given by Eq. 


Pr'ik) = Ph”{k) = pr*' 
driver nodes no for single layers with degree distribution P®' 
function of c and P(2). 


(60l. On the right: the doub le o f the density of 


(fc) = P°“‘(fc) = P{k) and P(fe) given by (601 is plotted as a 


E. Effect of degree correlations on controllability of multiplex networks 


In order to analyze the effect of degree correlations ^ on the controllability of multiplex networks, we correlate the 
degree of the replica nodes in the two layers of a duplex network formed by layer A and layer B. In particular we 
consider two cases: a duplex network in which only the low in-degree nodes (nodes of in-degree 0,1, 2 ) are correlated 
and a duplex networks in which the in-degrees of the replica nodes are correlated independently on their value. For 
each case, we define the joint in-degree distribution k^) between layers and the corresponding expression of 

the zero-temperature BP equations in the correlated ensemble of networks. 

In the first case we consider a joint in-degree distribution k^) given by 

r pStB k^P{k^) + {1- p)P{k^)P(k^), for k^ <2 
= J (I-p)P(/c^)P(fc^), for k^>2 fc^<2, 

\p^^P{k^) + (l-p)P{k^)P{k^), for k^>2 k^ > 2 

where C = 1 — J2k<2 where P{k) is a given normalized degree distribution. The distributions of the fields over 
the links of this ensemble of networks are given by 

= w'^6{h°‘-l)+w^6{h°‘+ l)+w^S{h°‘), 

= <d(/l“-l)-bu:2<5(^“ + l)+W3^(^“), (63) 


where a = A^B and where the probabilities wf, i 
Wi,W 2 ,'w^ that satisfy the equation Wi+tP^+w^ = 


iCg are normalized wf wf + u'f = 1 as well as the probabilities 
1. The zero-temperature BP (Max-Sum) equations (30) averaged 
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FIG. 11: Density of driver nodes e = no for a duplex network composed by two power-law networks with = 

= P(A:^’°“*) = P(fc) oc k~'^ and minimal degree m = 1 as a function of 7 (indicated with a solid red 
line). The minimum in/out degree 1 and the maximum in/out degree is given by the structural cutoff with N = 10"^. In the 
dashed blue line we display the double of the number of driver nodes e = 2 no for a single power-law network with the same 
in/out degree distributions P{k), indicating the fraction of driver nodes necessary to control separately the two layers. 


over this ensemble of networks can be expressed in terms of the probabilities {wf}i=i^ 2,3 and {wf }i=i, 2,3 as 


wi = p 


- (1 - + (Gi(u;2) - ^ - Go(l - tui)) 


(k) (k) 

-l-(l -p)Gi{w 2 ) [1 - Go(l - Wi)] 


W2 = p 




((.•) (t) 

P(l) 2P(2) 


(64a) 


(fc) 


tn n ^ 2P(2) ^ 

-(Gi(l - wi) - - ~ ~ ‘^0(^2)) 

-l-(l -p)[l- Gi(l - wi) -I- Gi(l - wi)Go {W2)] 


(64b) 
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FIG. 12: On the left: the density of driver nodes ud as a function of the parameters 7 and P(2) for duplex n etw orks of = 10® 
nodes with degree distributions = P^{k) = P^^(k) = P{k) and P{k) given by Eq. (611. On the right: 

double of the density of driver nodes no as a function of the paramete rs 7 and P{2) for single networks of = 10® nodes with 
degree distributions P^’^{k) — P°'^^{k) = P{k) and P{k) given by Eq. (611. 


where 

Go(z) = 

k 

G,{z) = Y. ( 65 ) 

k>3 


Finally the energy E is given by 

E = 2 {Go i'W2) - [1 - Go(l - til)]} + 2{k} [u)i(l - W 2 ) + wi(l - W 2 )] - 2(1 - p) {[1 - Go(l - tyi)][l - Go(w2)]} 


- 2p {P(l)u;i(l - W 2 ) + P(2)(l - (1 - u;i)2)(l - wl) + G(1 - Go(l - u;i))(l - Go(u; 2 ))} (66) 

In the second case we consider the joint degree distribution P™{k^, k^) given by 

P™(A:^, k^) = pSkB^k^Pik^) + (1 - p)P{k^)P[k’^), (67) 

where P{k) is a given normalized degree distribution. The distributions of the fields over the links of this ensemble 
of duplex networks are given by 

p„(/i“) = w“d(/i“-l)+ri;“<5(/i“ + l)+w“d(h“), 

V^(h°‘) = <d((i“-l)+ih“(5(/i“ + l)+u;“d(h“), (68) 

where a = A, B and where the probabilities wf,W 2 , are normalized wf + tcf + = 1 as well as the probabilities 

wf,11)2,11)^ that satisfy the equation + Wg = 1. We get the equations 

wi = p[Gi('u; 2 ) - (1 - wi)Gi(ui 2 (l - wi))] + (1-p)Gi(u; 2 ) [1 - Go(l - wi)] (69a) 

W2 = p[l-Gi{l-wi)+W2Gi{w2il-wi))] + {l-p)[l-Gi{l-wi) + Gi{l-wi)Go{w2)], (69b) 
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FIG. 13: Density of driver nodes no as a function of P{2) for a duplex network with PX'^k) — PX^^ik) = P^{k) = PjX'^{k) = 
P{k), P(k) given by Eq. (61 1 and 7 = 2.3. The fraction of driver nodes computed with the zero-temperature BP (Max-Sum) 
algorithm on a duplex network of = 10"^ nodes (averaged over 25 network realizations) is compared with the theoretical 
expectation for the density no in an ensemble of random duplex networks with the given degree distributions. 


where the generating functions Go{z) and Gi{z) are defined as 


Go(z)=^P(fc)z^ 

k 

= (70) 


The energy E in this ensemble is given by 

E = 2 {Go (W 2 ) - [1 - Go(l - thi)]} + 2{k) [■u:i(l - W 2 ) -I- wi(l - ti) 2 )] 

— 2(1 — p) {[1 — Go(l — tei)][l — Go(tC 2 )]} ~ 2p{l — Go(l — wi) — Gq{w 2 ) + Go('i« 2 (l ~ tt'i))} 


(71) 


The deg ree correlation of low in-degree nodes can modify the number of driver nodes nu found in duplex networks (see 


Fig. 14 for the case of a duplex network formed by Poisson layers with = c). 


Once the low in-degree nodes are correlated, correlating also the other in-degrees of the network does not change 
substantially the number of driver nodes ud as discussed in the main body of the paper. 
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FIG. 14: The density of the driver nodes no in a duplex network formed by two Poisson networks with = 

^ with correlated low in-degrees is plotted as function of c for different values of p. The result for 
the two separate layers is shown in black (dashed curve) while the situation for uncorrelated layers is shown in red. The value 
of p increases going from the red curve {p = 0 ) to the grey curve (p = 1 ). 
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